home *** CD-ROM | disk | FTP | other *** search
/ Info-Mac 3 / Info_Mac_1994-01.iso / Science / MathPad 2.1.5 / examples / curve fit < prev    next >
Text File  |  1993-06-21  |  1KB  |  39 lines
  1. ----------general least squares fit------
  2. -- finds coefficients a[i] for best fit of weighted sum of functions of x
  3.  
  4. fit(x) = sum(a[i,1]*f(x)[i],i,1,nparms)
  5.  
  6. f(x)[j]=x^(j-1)  -- polynomial a1+a2*x+a3*x^2...
  7.  
  8. P[i,j]=f(x[i])[j] dim[ndata,nparms]
  9. PtP:=multiply(transpose(P),P):
  10. Pty:=multiply(transpose(P),y):
  11. inv:=invert(PtP):
  12. a:=multiply(inv,Pty):
  13.  
  14. data=read(xydata)
  15. x[i]=data[i,1]; y[i]=data[i,2] dim[ndata]; 
  16.  
  17. ndata=count(data)
  18. nparms=3
  19.  
  20. a:{{5.582},{0.159},{0.006}}
  21.  
  22. plot data
  23. plot fit(X)
  24.  
  25. ---------- matrix operations ------------
  26. multiply(A,B)[i,j] = sum(A[i,k]*B[k,j],k,
  27.          1,count(B)) dim[count(A),count(B[1])]
  28. transpose(A)[i,j] = A[j,i] dim
  29.                         [count(A[1]),count(A)]
  30. invert(A) = adjoint(A)/det(A)
  31. adjoint(A) = transpose(cofactor(A))
  32. cofactor(A)[i,j] = (-1)^(i+j)*
  33.                    det(submatrix(A,i,j))
  34. submatrix(A,k,l)[i,j] = A[i,j] when i<k and j<l,
  35.             A[i+1,j] when i≥k and j<l,
  36.             A[i+1,j+1] when i≥k and j≥l,
  37.             A[i,  j+1] when i<k and j≥l dim
  38.                         [count(A)-1,count(A)-1]
  39.